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I. INTRODUCTION/BACKGROUND 



Two important parameters of interest in the field of aerodynamics 
of airfoils are lift and drag. Lift is defined as the force generated by 
an airfoil normal to the mean free stream direction of flow. Lift can be 
approximated by integrating the pressure distribution around an airfoil. 
Drag is the component of force parallel to the mean free stream direc- 
tion of flow. Drag on airfoils is usually defined to consist of two 
forms, skin friction drag and pressure drag where the total drag on 
the airfoil is usually defined as profile drag. Drag on airfoils is more 
difficult to calculate and predict than lift. 

Flow separation occurs when the fluid "particles" fail to follow the 
surface of the airfoil. Separation is sometimes confused with the 
phenomenon of stall. Stall can not occur without separation but separa- 
tion can occur without stall. Separation is usually a three dimensional 
effect and incorporating this effect into lift and drag calculations on 
airfoils which are two dimensional is a very difficult task. Prediction of 
drag on airfoils on a purely theoretical basis is difficult even when the 
flow has not separated since transition of the flow from laminar to 
turbulent significantly affects the drag. 

Laminar flow is distinguished from turbulent flow visually by the 
macroscopic motion of the fluid. Laminar flow appears smooth to the 
eye while turbulent motion is very erratic to the eye. The physical 
mechanisms that induce this motion are different. Generally speaking, 
most flows of interest start laminar (steady) and transition to turbu- 
lent (unsteady). The start and finish point of transition is very diffi- 
cult to predict. This transition point or range is usually empirically 
determined . 

The unsteady nonlinear three dimensional Navier Stokes equations 
define the macroscopic motion of a Newtonian fluid. A Newtonian fluid 
is a fluid in which the stresses are related to the rates of strain of the 
fluid particle and the coefficient of viscosity. Air is assumed to be 
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governed by these equations. However, until recently there was no 
practical method of solution to these equations. Due to this difficulty 
the Navier Stokes equations were reduced to two general sets of 
equations : 

1. I nviscid/potential 

2. Viscous/boundary layer. 

The inviscid equations are further simplified by assuming steady 
flow, no external body forces, and no energy transfer to or from the 
fluid. With the additional assumption of uniform approach flow the flow 
is irrotational and hence permits the introduction of a velocity poten- 
tial. The inviscid/potential flow solution is of limited value because it 
predicts no viscous drag (D'Alembert's paradox). However, in combina- 
tion with the Kutta-Joukowski condition the inviscid solution predicts 
the lift of airfoils with sharp trailing edges at small angles of attacks 
very well. The condition originally proposed by Kutta in 1902 and 
independently by Joukowski in 1906 [Ref. 1: pp. 390-395] asserts only 
that the upper and lower surface flows meet smoothly at the trailing 
edge. This condition is sometimes referred to as a zero load condition at 
the trailing edge, meaning there is a pressure equalization at the 
trailing edge. 

The Kutta condition is a practical solution to the problem that the 
inviscid flow solution is singular when the flow negotiates a sharp 
corner (that is the velocity goes to infinity) for an arbitrary circula- 
tion. Placement of a rear stagnation point eliminates this problem and 
introduces a unique solution to the inviscid equations. Experimental 
evidence supports the theory of circulation and the Kutta condition on 
airfoils with relatively sharp trailing edges at small angles of attack. 
However, in certain flow regimes the inviscid/potential flow solution 
with the Kutta condition fails to adequately explain real world flow 
phenomena, such as flows around wings at high angles of attack and 
practical turbomachinery blades in cascades. Here the effects of 
viscosity can no longer be ignored. 

With the effects of viscosity included, the unsteady Navier Stokes 
equations are difficult to solve. Simplification is accomplished by 
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discarding higher order terms. These simplified Navier Stokes equa- 
tions are now termed boundary layer equations. The boundary layer 
concept , first introduced by Prandtl in 1904, defines a very thin 
region next to the surface of an object where all significant velocity 
gradients normal to the streamwise direction occur and hence viscous 
forces are generated. Outside this region all viscous forces are 

assumed negligible and the inviscid/potential solution is assumed to 
adequately model the flow. The boundary layer equations can predict 
the drag on airfoils since most drag is due to skin friction caused by 
the viscous forces in attached flow. The problem is made more difficult 
by the need to compute both laminar and turbulent boundary layers and 
to determine the transition point from laminar and turbulent flow. The 
reader is also reminded that the boundary layer concept is applicable 
only at sufficiently high Reynolds numbers. For a more detailed 
discussion of the boundary layer concept see Schlichting [Ref. 2] and 
Cebeci and Bradshaw [Ref. 3]. 

Turbulent flows are difficult to predict and model. Alleviation of 
the unsteady nature of turbulent flows is generally accomplished by 
expressing the pertinent variables of the unsteady three dimensional 
Navier Stokes equations as a time averaged quantity and a fluctuating 
quantity and then time averaging the equations term by term. These 
time averaged Navier Stokes equations produce nine additional unknown 
entities that can be correlated to a stress tensor. These extra turbu- 
lent unknowns are termed Reynolds stresses. However, all nine 

Reynolds stresses are not considered significant. In fact, only one 
Reynolds stress is usually retained in a two dimensional analysis, which 
is the case for the Cebeci program used in this thesis. This Reynolds 
stress is the shear stress that is in the normal plane to the surface in 
the streamwise direction. It is considered the turbulent counterpart to 

the laminar shear stress. This particular Reynolds stress is sometimes 
modeled as the product of an empirical "eddy viscosity" and the local 
velocity gradient. 

Computational fluid dynamics is becoming a more important field as 
memory capacities grow bigger, computers get quicker, and numerical 
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algorithms get better. Accurate prediction of flows minimizes the 

amount of wind tunnel testing and helps the designer to make better 
designs. Two "state-of-the-art" programs for predicting flows and 
boundary layers in a cascade the Q3DFLOW-81 and a CEBCAS are exam- 
ined. The Q3DFLOW-81 is a turbomachinery program developed by Dr. 
Hirsch and associates of the Free University of Brussels [Ref. 4] and 
[Ref. 5]. The blade-to-blade portion is used to calculate inviscid 
incompressible flows in a cascade. Inviscid results are presented for 
several airfoils/blades and configurations. Experimental data, analytical 
solutions or other computational flow data are compared where appli- 
cable. The CEBCAS is a standard Cebeci two-dimensional, incompres- 
sible and constant viscosity boundary layer program modified to give 
graphical output for boundary layer and skin friction distribution and 
calculate losses and drag coefficients for airfoils/blades in cascade. 
The CEBPLT is also a standard Cebeci two-dimensional, incompressible 
and constant viscosity boundary layer program modified to give graph- 
ical output for boundary layer and skin friction distribution and calcu- 
late drag coefficients for isolated airfoils. CEBPLT and CEBCAS are 
used to calculate boundary layers on several isolated airfoils and 
airfoil/blades in a cascade. Experimental data is correlated when avail- 
able and reasonable to do so. Additionally, the theory and operation of 
the Cebeci boundary layer programs (Chapter III. and Appendix A 
respectively) are presented. 
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II. I N VI SC I D FLOW COMPUTATIONS IN A CASCADE 



Physical flow is three dimensional and unsteady, especially in an 
axial flow compressor. Simplification of this complex flow in a 
compressor is usually approached by splitting it into two two- 

dimensional flows: 

1 . The meridional flow 

2. The cascade flow (sometimes referred to as the blade-to-blade) 

Solution to the flow in the cascade plane is input for the meridional 
"through flow" solution, which makes the accurate prediction of flow in 
the cascade plane very important. Gostelow [Ref. 6: p. 86] has attrib- 
uted systematic advances in blade design to the ability to calculate 

flows in a cascade. 

Inviscid flow computations in a cascade are indeterminate for prac- 
tical turbomachine blades. The outlet condition, usually the outlet 
angle, is specified and at present is found experimentally. For isolated 
airfoils the Kutta condition gives meaningful results at small angles of 
attack for practical shapes (that is airfoil shapes that have relatively 
sharp trailing edges). However, most practical turbomachinery blades 
have trailing edges that are rounded so as to avoid stresses and manu- 
facturing problems. Frequently, laminar flow will separate near the 
leading edge and reattach turbulently downstream. This further 
complicates the modelling process of the flow. Using only an inviscid 
model it is very difficult to predict the flow. Miller and Serovy 
[Ref. 7] examined five different trailing edge hypotheses to determine 
the deviation angle for cascades using the inviscid model. None of the 
five trailing edge conditions tested were considered acceptable as a 
general rule such as the Kutta condition is for sharp trailing edges. 
This is understandable since the outlet angle is influenced by viscous 
effects and on a round trailing edge the determination of a rear stagna- 
tion point is not as easy as on a sharp trailing edge. Gostelow 

[Ref. 8] reported in a particular example solved by exact potential flow 
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theory that a change in rear stagnation point location of as little as .3% 
of the chord would change the outlet angle by as much as 10 degrees. 
Also, the precise location of the rear stagnation point at some incidence 
is not necessarily the "correct" location for some other incidence. 

A. Q3DFLOW-81 (BLADE-TO-BLADE MODULE) 

The Q3DFLOW-81 is a turbomachinery program for calculating flows 
in an axial flow compressor or turbine. This program uses a finite 
element method. The blade-to-blade module was used to solve the 
inviscid equations for flow in a cascade. Substantial runs were made 
with the Q3DFIOW-81 program for various compressor blades (including 
even sharp trailing edge airfoils) in a cascade. In order to compare 
with experimental data each configuration was run twice, one time with 
outlet angle specified by the experimental angle and another with the 
Kutta condition option. 

The major problem encountered in running the Q3DFLOW was in 
generating an acceptable finite element mesh. The finite element mesh 
generation appeared to be very sensitive to the placement of coordinates 
used to define the shape of the airfoil/blade, especially around the 
leading edges and trailing edges. Observations of running Q3DFLOW 
(blade-to-blade) indicate a need for a close one-to-one correspondence 
between lower and upper coordinates of the surface of the airfoil, say 
within one to three percent of chord (i.e. the streamwise direction) for 
the major portion of the airfoil and closer for leading and trailing edge 
if rounded. The more curvature in the airfoil geometry the smaller the 
spacing between adjacent coordinates needed to define the geometry to 
ensure that flow computations are for the right shape. However, too 
many points to define the leading edges or rounded trailing edges will 
generate erroneous mesh elements. ' A general rule of thumb is to use 
the minimum points needed for the eye/brain to generate the correct 
airfoil shape by mentally connecting coordinates. 

Four parameters ILWIN, IVPI N, I LWOUT, and IVPOUT of CARD NR 
B-03 of users guide [Ref. 9] are very useful in solving finite element 
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mesh problems. These parameter modified the mesh around the leading 
and trailing edges by forming a triangle element instead of the normal 
rectangular shaped elements ( Figure 2.1). Two other parameters that 
are sometimes useful are NSTUP of CARD NR B-01 (number of stations 
in the upstream) for solving leading edge mesh problems and NSTDWN 
of CARD NR B-01 (number of stations in the downstream) for trailing 
edge problems (Figure 2.1). A "last ditch maneuver' is to input 
different coordinates to define the airfoil/blade shape. 
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III. THEORY OF CEBECI BOUNDARY LAYER PROGRAM 



The governing equations (Eqs. 1, 2 and 3) of the Cebeci two 

dimensional boundary layer program are obtained from an engineering 
simplification of the two dimensional Navier Stokes equations by 
neglecting small terms due to the known empirical fact that the 
displacement thickness is very small compared to the characteristic 
length of the airfoil. These equations are then assumed to be appli- 
cable to turbulent flows by including only the Reynolds shear stress 
and neglecting the two Reynolds normal stresses that are obtained by 
time averaging the two dimensional Navier Stokes equations. Reynolds 
stresses are essentially empirically obtained and can be expressed in 
many different forms, depending on the theory used. 




Eq . 1 
Eq. 2 
Eq . 3 
Eq . 4 



The eddy viscosity formulation is used in the Cebeci program to 
model the Reynolds shear stress. This concept was first introduced by 
Boussinesq in 1877 [Ref. 2: p. 578] in which the turbulent shear 

stresses are assumed to be similar in nature to the laminar shear 
stresses and are therefore proportional to the velocity gradient (Eq. 
4). Cebeci and Smith [Ref. 10: pp. 255-257] have formulated a two- 
range empirical eddy viscosity model. In the inner region the formula 
is similar to the Prandtl mixing length formula modified by a Van Driest 
damping parameter and intermittency factor that accounts for the tran- 
sitional region. The outer region is modelled as proportional to an 
empirical constant, the velocity at the edge, the displacement thickness, 
and also to the intermittency factor. 
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Assuming constant density and viscosity and using definitions as 
stated in Eqs. 5 and 6, Eq . 1 becomes Eq. 7. For laminar flows b is 
equal to one. Also b may be redefined so that Eq. 7 is useful for 
axisymmetric laminar flows (Mangier Transformation) with small tran- 
sverse curvature effect [Ref. 3: p. 48-50]. 
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b=l+e* 
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Eq . 6 
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Since the rate of change of local static pressure in the y direction 
is simplified to equal zero (Eq. 3) and the velocity at the edge of the 
boundary layer (that is when y is equal to the displacement thickness) 
is a function of x only, then the Bernoulli equation (Eq. 8) may be 
used to determine the rate of change of local static pressure in the x 
direction, as the value in the free stream and surface are considered 
equal. Differentiating Eq. 8 with respect to x gives Eq. 9. (Note that 
Eq. 9 is now a differential vice partial since the local static pressure is 
a function of x only.) The local static pressure is now specified as a 
pressure coefficient, velocity ratio, etc. from the inviscid/potential 
solution as a boundary condition for the boundary layer equations. 
With Eq. 9 ,the simplified momentum equation in the x direction (Eq. 7) 
and the continuity equation (Eq. 2) reduce to a system of two partial 
differential equations with two unknowns u and v for laminar flows. 



P + 2P u »"const 
do _ du. 



Eq. 8 
Eq . 9 



The Falkner-Skan transformation (Eq. 10) defines a dimensionless 
coordinate ti and Eq. 11 defines a stream function which consist of a 
dimensionless stream function f in x and n. A necessary and sufficient 
condition for the existence of a stream function is that continuity (Eq. 
2) be satisfied and hence Eqs. 12 and 13 must be satisfied. This 
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transformation along with the definition of a stream function and use of 
the chain rule for differentiating reduces Eqs. 7 and 2 into one partial 
differential equation with one unknown stream function in ri and x for 
laminar flows with no known analytical solution . For similar flows this 
transformation changes a system of two partial differential equations 
into one ordinary differential equation with an unknown stream function 
as a function of the independent variable n . 



<\J 

s 

-X 

in 


Eq. 10 


^=(u e ux) I/2 f (x,^) 


Eq. 11 


“-$5 


Eq. 12 




Eq. 13 



Similar flows are laminar flows that satisfy Eq. 14. An example of 
a similar flow is the flow over the flat plate. In this example the 
pressure-gradient parameter, m is equal to zero. The definition of m 
(Eq. 15) may be derived by differentiating Eq. 14 by x substituting 
Eq. 14 back into this differentiated equation, and then solving for m. 
The concept of similarity is not as important now as it was before the 
advent of the computer. It allowed for practical numerical solutions of 
certain flows even though very few of the flows can model any real 
world phenomenon. For further understanding of the concept of simi- 
larity see References 2 and 3. 



u 




m 

x 



m - 



X— c^u. 
u e 3x 



Eq . 14 
Eq . 15 



Derivation of Eq. 22 is a simple matter of putting u , v and their 
partial derivatives in terms of the dimensionless stream function f by 
using Eqs. 10-13 and the chain rule. (For example Eq. 16 gives Eq. 
17, where prime denotes partial differential with respect to ri). Eqs. 
17-21 show the results of this. Substituting Eqs. 17-21 into Eq. 7, 
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rearranging , expanding, and multiplying the equation by x and 
dividing by u e squared, and then substituting in the pressure gradient 
parameter (Eq. 15), Eq. 22 may be obtained for non-similar laminar 
flows (b=1). Now b is a function of n and must be 
to account for turbulent flows. 
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This third order nonlinear partial differential equation (Eq. 23) is 
written as a system of three first order partial differential equations 
(Eq. 24-26) by defining two new variables script u and script v func- 
tions in n (Eq. 24 and 25 respectively). Script u is a nondimensional 
x-component velocity but script v has no direct correspondence to the 
y-component velocity. The three unknowns of this third order system 
are script u, script v, and f. Recalling Eq. 17 and the fact that the 
velocity is zero at a solid surface of an airfoil ( n =0) , the stream func- 
tion f and its partial derivative with respect to n is equal to zero. 
Hence two of three boundary conditions needed are Eq. 27 and 28. 
The third boundary condition (Eq. 29) comes from the velocity at the 
edge (n = displacement thickness) being modelled equal to the velocity at 
T| equal to infinity for a given x location. 
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order in n on the net rectangle shown in Figure 3.1 and averages at 
the midpoints of the net. The central finite difference of the first 
derivative formula (Eq. 32) used in Cebeci's program can easily be 
derived by use of a Taylor series expansion of a point x + h and a 
point x-h (Eq. 30 and 31 respectively , see Figure 3.2). Eq. 32 is 
obtained by subtracting Eq.31 from Eq.30 and then solving for the first 
derivative. 
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Figure 3.1 Boundary Layer Grid 
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f(x+h)-f(x)+hf (x)+hlgher order terms 

f(x-h)-f (x)-hf (x)+htgher order terms 

r '( x )-L (x* h ?; h r( x -h) l0 (h2) 

The 0(h 2 ) in Eq. 32 represents the sum of the remaining terms of order 
h 2 in the series. Eq. 29 gives an exact value of the derivative if all 
terms in the infinite series are retained or if h goes to zero. Since 
Eq. 29 is truncated by discarding all terms of order h 2 this finite 

difference formula is of truncation error of order h 2 . The spacing 

between grid points (h) is usually relatively small and the error 
incurred by this approximation is usually insignificant for most engi- 
neering purposes. Central finite differences of the first derivatives are 
applied at the mid points of the net. Twice the hx and hy distance are 
defined as kn and hj (that is the x distance between i-1 and i grid 
location and the ri distance between j-1 and j grid location respec- 
tively, see Figure 3.1). Rewriting Eq. 32 in terms of subscript 

sequence grid location numbers i and j of the spanwise and normal 
direction respectively at the midpoints and accounting for the fact that 
all derivatives are partial gives Eq. 33 and Eq. 34 respectively. For 
example, Eq. 33 is the partial derivative of f in the spanwise direction 
for any constant normal grid location (j). Using Eq. 33 to approximate 
Eq. 24 and 25 averaging script u and script v at the midpoint in the j 
direction , and rearranging, Eqs. 35 and 36 are obtained. In a similar 
manner and with use of the central difference approximation for the 
partial derivative in the spanwise direction Eq. 26 is transformed from a 



partial differential 


equation to an algebraic 


equation and will 


not be 


presented 


here due 


to its lengthy derivation 


and additional definitions 


introduced 


that are 


unnecessary but convenient [Ref. 3: pp. 


216 and 


217]. 


f' - f 


\ 

'-'•'J 


Eq . 33 






/2, j 


k 
_ * 






f' - * 


wj ' ifj-l 


Eq . 34 






j-l/2 


h i 





Eq . 30 
Eq . 31 
Eq. 32 
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Eq . 36 



Eq . 35 







X-h x X + 'n 



Figure 3.2 Points Defined for Taylor Series Expansion 
If all values of f, script u, and script v are known at some initial 



35 , 36 and Eq. (7.2.10c) of Reference 3 are a nonlinear system (that 
is the system is nonlinear since the last equation is nonlinear) of 3J + 3 
equations with 3J + 3 unknowns where j=0, 1 ,2,3. . . , J . Cebeci uses 
Newton's method to solve this system of nonlinear equations. Cebeci 
departs slightly from the strict application of Newton's method in 
treating the b terms of Eq. (7.2.10c) of Reference 3 to avoid unneces- 
sary complication [Ref. 11]. The solution to this now linearized set of 
equations is accomplished by a very efficient block elimination method 
which is possible due to the unique form that the system of equation 
can be put in. This method is fairly elaborate and will not be 
presented here. 

The basic idea of Newton's method can be understood from a little 
geometry. Suppose we wish to solve for the roots of this nth order 
polynomial in the independent variable x where n is greater than 1 and 
hence nonlinear. Say this polynomial f(x) is defined as in Figure 3.3 . 



downstream spanwise grid location, say i-1 for all values of j then Eq. 
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Now suppose we guess a solution to one of the roots of this polynomial 
say xl . From geometry, the tangent of 0 is f(x1)/x1-x2 which is also 
the slope at xl that is f'(xl). Solving for x2 and writing in terms of 
sequence iteration numbers for successive approximations we obtain Eq. 
37. An exact value to this root can not be given with only a finite 
number of iterations so an acceptable convergence limit is usually 
defined between successive iterations. Cebeci considers the wall shear 
parameter that is f" at ti equal to zero (also known as script v) as the 
greatest source of errors in boundary layer equations. The conver- 
gence criterion is .00001 between successive iterations of the wall shear 
parameter (f"(0)) for laminar flows. For turbulent flows the conver- 
gence criterion is expressed as a percentage. There are many problems 
that can be encountered by this technique. For instance if the initial 
guess is "poor" this technique may converge to a root not sought or 
may not even converge at all. However, due to the nature of the 
boundary layer equations and the manner in which Cebeci applies it, 
Newton’s method is "unconditionally stable" For a more detailed discus- 
sion of the theory of Cebeci boundary layer program see Reference 3. 










xZ 



xi 



Figure 3.3 Newton's Method 




Eq. 37 
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IV. RESULTS OF Q3DFLOW INVISCID BLADE-TO-BLADE 

COMPUTATIONS 



Several blade shapes and configurations were run using the blade 
to blade portion of the Q3DFLOW-81 program. Several different graph- 
ical outputs are available from Q3DFLOW, but only those that could be 
compared with other available information are presented here. Results 
are presented in graphical output usually as pressure coefficients. 

A no camber, sharp trailing edge airfoil, the NACA 65-010 (see 
Appendix B for shape and coordinates used to input Q3DFLOW) was run 
at an inlet (3 = 30°, a = 1.5, and several a's. Figure 4.1 shows the 
results of the Q3DFLOW-81 blade to blade (solid line) plotted as a pres- 
sure coefficient vs. chord. The circle and triangle points depict the 
experimental pressure coefficients obtained from NACA TN 3916 
[Ref. 12]. Figure 4.1 shows the results specifying the experimental 
outlet 3 as given in NACA TN 3916. Figure 4.2 shows results using 
the Kutta condition. 

A slightly cambered, rounded trailing edge blade, the NACA 65-(4Aj, 
I 0b )10 (see Appendix C for shape and coordinates used to input 
Q3DFLOW) was run at an inlet S =30°, a =1.0 and several a’s. Figure 

4.3 shows the results of Q3DFLOW-81 blade to blade (solid line) as 

pressure coefficient vs. chord when specifying the experimental outlet & 
as given in NACA TN 3817 [Ref. 13]. Figure 4.4 shows the results 

when using the Kutta condition. Experimental pressure coefficients 
(circles and triangles) obtained from NACA TN 3817 are also plotted. 

A slightly cambered, sharp trailing edge airfoil, the NACA 8410 

(see Appendix D for shape and coordinates used to input Q3DFLOW) 

was run at an inlet 3 = 44.4° and a = 1.33. Figure 4.5 shows the 
results of Q3DFLOW-81 (solid line) using the Kutta condition, as velocity 
ratios referenced to the approach velocity vs. chord. The velocity 
ratios, as calculated by the Schlichting-Method and given by 
AGARD-AG-220 [Ref. 14: pp. 317-319], are referenced to the approach 



26 



REPRODUCED AT GOVERNMENT EXPENSE 







0^0 CN0B0 





Figure 4.1 Q3DFLOW Results for NACA 65-010 
(Outlet 6 Given from Experiment) 
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Figure 4.2 Q3DFLOW Results for NACA 65-010 
(Kutta Condition) 
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Figure 4.3 Q3DFLOW Results for NACA 65- (4A.L. ) 10 
(Outlet 3 Given from Experiment) 



29 



PRESSURE COEFFICIENT PRESSURE COEFFICIENT PRESSURE COEFFICIENT 



REPRODUCED AT GOVERNMENT EXPENS 









Figure 4.4 Q3DFLOW Results for NACA 65-(4 A,Iol)10 

(Kutta Condition) *■ 5b 
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velocity vice the mean velocity as originally given and are presented as 
circles and triangles in Figure 4.5 . 

A cusped blade (see Appendix E for shape and coordinates used 
to input Q3DFL0W) from Gostelow [Ref. 6: pp. 109,110,123,124] was 
run at an inlet 0 = 37.5° and o =1.01. Figure 4.6 shows results as 
pressure coefficient vs. per cent chord using the Kutta condition. The 
circles and triangles on Figure 4.6 show the results as analytically 
determined by Gostelow. 

An advanced controlled-diff usion fan blade (see Appendix F for 
shape) was run at an inlet 3 = 68°, X - 58.25° ,o = 1.34 and an axial 
velocity/density ratio ( AVDR) = 1 .443 . Figure 4.7 shows results of 
Q3DFLOW-81 as pressure coefficient vs. chord (solid line) using the 
experimental outlet angle as measured in the cascade wind tunnel at the 
Turbopropulsion Laboratory of the Naval Postgraduate School. Figure 
4.8 shows results using the Kutta condition. Experimental pressure 
coefficients as measured in the cascade wind tunnel at the 

Turbopropulsion Laboratory of the Naval Postgraduate School are 
plotted as circles and triangles. 

Table I shows outlet 3 as calculated from Q3DFLOW using the Kutta 
condition option for all blades run. The Kutta condition option of 
Q3DFLOW uses an initial outlet 3 which was given as the experimental 
outlet 3 or other analytical outlet 3. Also shown is outlet 3 given from 
experiment or another available method. 
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Figure 4.5 Q3DFLOW Results for NACA 8410 
(Kutta Condition) 
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Figure 4.6 Q3DFLOW Results for Cusped Blade 
(Kutta Condition) 
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Figure 4.7 Q3DFLOW Results for Fan Blade 
(Outlet 3 2 Given from Experiment) 
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Figure 4.8 Q3DFLOW Results for Fan Blade 
(Kutta Condition) 
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TABLE I 

Outlet Angle Comparison 



NACA 


65-010 


INLET 6=30° 


a= 1.5 






a 


Q3DFLOW 
EXPERIMENTAL KUTTA 


% DIFF 




-3.° 


33.9° 


34.77° 


2.6 




2.° 


28.7° 


29.48° 


2.7 




4.° 


26.6° 


27.37° 


2.9 




6.° 


24.6° 


25.27° 


2.7 




9.° 


21.7° 


21.77° 


.3 




15.° 


16.2° 


16.18° 


.1 


NACA 


65- (4A 


1 )10 INLET 6=30° o=1.0 






a 


EXPERIMENTAL 


03 D FLOW 
KUTTA 


% DIFF 




-2.4° 


30.2° 


29.18° 


3.4 




3.4° 


25.9° 


23.43° 


9.5 




5.4° 


21.9° 


20.73° 


5.3 




8.4° 


18.9° 


17.93° 


5.1 




11.4° 


16.8° 


14.83° 


11.7 




17.4° 


11.1° 


9.18° 


17.3 


NACA 


8410 


INLET 6=30° c 


>=1 .0 








SCHLICHTING 


Q3DFLOW 

KUTTA 


% DIFF 






21 .3° 


20.26° 


4.9 


FAN 


INLET 6=30° o=1.0 


AVDR=1 


.43 






Q3DFLOW 
EXPERIMENTAL KUTTA 


% DIFF 






48.62° 


41.66° 


14.3 


CUSPED BLADE INLET 6=30 


o=1.0 








GOSTELOW 


Q3DFLOW 

KUTTA 


% DIFF 






30.02° 


30.02° 


0.0 
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V. RESULTS OF CEBECI BOUNDARY LAYER PROGRAM 



Several isolated airfoils were run on CEBPLT. Several airfoils and 
turbomachinery blades were also run using pressure or velocity distri- 
butions from inviscid flow in a cascade (Q3DFLOW-81 blade to blade) as 
input for CEBCAS. Results from CEBPLT and CEBCAS are presented 
as graphical output displaying displacement, momentum and skin friction 
distributions. The velocity distribution input is also displayed. 
Transition was computed by Michel's method as given by Cebeci and 
Smith [Ref. 10: pp. 332-333] or the point of laminar separation if this 
occurs first. Drag coefficients on isolated airfoils were calculated using 
a Squire- Young formula [Ref. 2: pp. 764-765]. Losses were calculated 
for airfoils/blades in a cascade as given by Roudebush and Lieblein 
[Ref. 15: p. 24]. The drag coefficient for cascades was calculated by 
multiplying the losses by the cosine 0 at infinity (i.e. the 0 associated 
with the mean free velocity) and dividing by a. 

NACA 0012 was run using velocity ratios obtained from an inviscid 
program of the Douglas Aircraft Company at a = 0°. NACA 0009 and 
NACA 0006 were run using inviscid velocity ratios from Abbott and von 
Doenhoff [Ref. 16] also at a = 0°. NACA 4412 was run using the 

experimental pressure coefficients from Pinkerton [Ref. 4] again at a = 

0° . Figures 5. 1-5.5 show the results of these airfoils as displacement 
thickness, momentum thickness, and skin friction distributions vs. 

x/c, fraction of chord. The profile drag coefficient , C^ , and skin 
friction drag coefficient, C + Referenced to approach velocity and 
velocity at the edge, u e are for one surface only and are given on the 
skin friction graph. The total profile drag for symmetrical airfoils is 
twice the C value given on skin friction graph at a = 0°. Table II 
shows the comparison between calculated profile drag and measured 
profile drag as given in [Ref. 16]. 

The NACA 8410 airfoil was run on CEBCAS using the velocity 

distribution as given by the Schlichting method in [Ref. 14: pp. 
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Figure 5.1 NACA 0012 Boundary Layer Results 
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Figure 5.2 NACA 0009 Boundary Layer Results 
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Figure 5.3 NACA 0006 Boundary Layer Results 
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Figure 5.4 NACA 4412 Boundary Layer Results (Suction Side) 
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Figure 5.5 NACA 4412 Boundary Layer Results (Pressure Side) 
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317-318]. after referencing to the approach velocity as input at a 
Reynolds number of 500,000. Figures 5.6 and 5.7 show results as 
thickness and skin friction distributions with transition specified as 
close to the leading edge as possible without incurring a program error. 
Also shown are the velocity distribution inputs. The displacement and 
the momentum thickness as calculated with the help of Truckenbrodt’s 
method for a fully turbulent boundary layer as given in [Ref. 14] are 
plotted as squares and circles, respectively, on the thickness distribu- 
tion graphs. 

The NACA 65-010 and the NACA 65-(4A a lg b )10 were run on 
CEBCAS using the inviscid pressure distribution from Q3DFLOW as 
input at Reynolds number of 245,000 and 444,000, respectively. Figure 
5.8 and 5.9 show results as thickness and skin friction distributions for 
the NACA 65-010 at an a = 6°. The drag coefficient for the cascade 
referencing to approach velocity as calculated by CEBCAS for the 
NACA 65-010 at an a = 6° was .0138 as compared to .013 as measured 
in NACA TN 3916. Table III shows where separation occurs according 
to CEBCAS in x/c. A "NO" indicates that separation did not occur 
over forward 95% of chord. Lieblein's equation used for losses does 
not account for separated flows or blades with round trailing edges and 
for that reason no other comparisons were made with experimental data. 
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Figure 5.6 NACA 8410 Boundary Layer Results (Suction Side) 
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Figure 5.7 NACA 8410 Boundary Layer Results (Pressure Side) 
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Figure 5.8 NACA 65-010 Boundary Layer Results (Suction Side) 
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Figure 5.9 NACA 65-010 Boundary Layer Results (Pressure Side) 
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TABLE II 

Comparison of Calculated and Measured Drag Coefficients 



AIRFOIL 


a 


REYNOLDS 


PROFILE 

EXPERIMENT 


DRAG 

CALCULATED 


NACA 0012 


0° 


6 MILLION 


. 00600 


. 00562 


NACA 0009 


0° 


6 MILLION 


.00575 


.00614 


NACA 0006 


0° 


6 MILLION 


.00510 


.00610 


NACA 4412 


0° 


3 MILLION 


. 00584 


. 00693 



TABLE III 

Separation Points from CEBCAS 

NACA 65-010 (Inlet B=30° o=1.5) 

Approximate x/c 



a 


Upper 


Lower 


-3° 


.7 


.025 


2° 


.65 


.025 


4° 


.65 


.025 


6° 


NO 


NO 


9° 


.05 


.6 


15° 


.025 


.65 


NACA 65(4A 


1 )10 


(Inlet 3=30‘ 




Approximate x/c 


a 


Upper 


Lower 


-2.4° 


.85 


NO 


3.4° 


NO 


.01 


5.4° 


NO 


.01 


8.4° 


NO 


.6 


11.4° 


NO 


NO 


17.4° 


.025 


.65 
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VI. DISCUSSION/CONCLUSIONS 



A. DISCUSSION 

The experimental pressure coefficients as found by the program 
Q3DFLOW agreed well with experimental data when using the Kutta 
condition or inputting the experimental outlet 3 with the exception of 
the fan blade. Q3DFLOW also gave good results when compared with 
the Schlichting and Gostelow analysis. Disagreement between experi- 
ment and Q3DFLOW for the fan blade case may have been due to the 
fact that the incidence angle run was very close to complete stall (i.e. 
wind tunnel stall) as evidenced by subsequent runs in the cascade 
wind tunnel at the slightly higher incidences. The axial velocity/ 
density ratio as measured in the wind tunnel was very high in this case 
and indicates lack of two-dimensionality. Q3DFLOW accounts for this as 
a blockage factor through the cascade but is only an linear approxima- 
tion. These two factors may have contributed to the disagreement 
between the Q3DFLOW results and experimental data as given for the 
fan blade. Q3DFLOW may give unrealistic pressure coefficients also . 
One encountered was as high as -34 . The coordinates used to define 
the leading and trailing edges may have been the cause since a judicial 
choice of the coordinates generally alleviated the unrealistic pressure 
coefficients . 

The discontinuous pressure coefficient distribution at the trailing 
edge on NACA 65-(4A l l ftVj )10 (Figures 4.3 and 4.4) is probably due to 
the discontinuous geometry at the trailing edge that was input to 
Q3DFLOW as a "last ditch maneuver" to generate an acceptable finite 
mesh. Initial attempts to generate an acceptable finite element mesh 
using coordinates where the trailing edge was sufficiently defined were 
unsuccessful. (Figure 6.1, coordinates were obtained by blending in the 
forward 95 % of blade with a 1 % chord radius trailing edge located as 
given in NACA TN 3817.) NACA TN 3817 indicates that the NACA 
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65- (4A ^ l 0b ) 10 may have possibly separated on the concave surface 
(pressure side) at a's of 3.4° and 5.4° for inlet (5 = 30° and a = 1. 
Figures 4.3 and 4.4 did not show any significant difference between 
measured and calculated pressure coefficients. 




Figure 6.1 NACA 65- 
Twith Round Trailing 

Experimental data on boundary layers in cascades is sparse at 
best which makes it difficult to evaluate the merits of the Cebeci 
boundary layer program. A further difficulty is that most cases run on 
the Cebeci boundary layer program using Q3DFLOW inviscid incompres- 
sible pressure distribution as input indicated separation. High pres- 
sure gradients and low Reynold numbers on the leading edge may make 
the boundary layer equations not valid. However, the calculated drag 
coefficient for the one case that had a sharp trailing edge and had 
attached flow on both sides, the NACA 65-010 at a = 6°, did agree 
well with the measured value. 

B. CONCLUSION 

Results of the program Q3DFLOW indicate that reasonably correct 
pressure distributions may be calculated when inputting the leading 
edge and trailing edge coordinates of a practical compressor blade with 



(4AJ flb )10 Sha 
Edge 8 Defined) 



Shape 
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care. The usefulness of the Cebeci program in calculating boundary 
layers in cascades still remains a question. Use of these two programs 
in conjunction with each other will probably not give accurate aerody- 
namic parameters such as lift and drag on practical turbomachinery 
blades in general due to separation being calculated by the Cebeci 
boundary layer program. 
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APPENDIX A. 

OPERATION OF CEBECI BOUNDARY LAYER PROGRAM 



The boundary layer program as sent from Cebeci in September, 1984 
calculates total skin friction drag, total profile drag using Squire- 
Young formula, and transition based on an empirical formula given by 
Cebeci and Smith [Ref. 10] in addition to those parameters presented in 
Cebeci and Bradshaw [Ref. 3]. Additionally, the inputs have been 
modified. The current version of this program used for this thesis has 
been slightly modified again to give graphical output and to handle 
numerous data files at once. Another slight modification can calculate 
losses for cascades using an equation 44 of Reference 15. 

A. INPUTS 

All inputs are formatted by fortran format. A card will refer to an 
input line that is a maximum of 80 columns. Each card will be defined 
by columns, format, variable and explanatory comment. For example: 

col 1-5 15 variable explanatory commentary 

col 1-5 indicates the field assigned to hold 

the value of the variable 

15 indicates the format used 

variable indicates a fortran name for the 

variable 

explanatory comment gives information on the variables 



Integer format statements are right justified that is in the preceding 
example if you want the variable to have a integer value of 1 you would 
put 1 in column 5. If you put the 1 in column 3 the value would be 
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read as the integer 100. The other two fortran format statements used 
in this program are F10.0 and 20A4. The F10.0 refers to a floating- 
point variable. A real number must be specified within the ten position 
field using a decimal point to separate the decimal point from the 
mantissa. Any amount of numbers right of the decimal point may be 
specified as input as long as it fits in the field. The 20A4 refers to a 
Hollerith variable. 80 positions are available for 80 characters. This 
format is used to read in the title for the output and the disspla 
routine (graphics). The last two character of the title must be S' so 
the disspla routine may read your title. 80 characters are available to 
specify the title but twenty to thirty are only recommended for the 
current character height used in the program. This may easily be 
modified for those familiar with DISSPLA ( a graphics software 
package) . 



1 . Inputs for C^BPLTj 



The following cards give the format to the modified Cebeci program 
(CEBPLT) that gives graphical output but does not calculate losses for 
cascades . 

$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$s$$$s$$$$$s$$ 

Run card Format(l5) 

********** 

col 1-5 15 N the number of runs (data files) 

$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$S$$$$$SSS$SSS$SSSS$$$SS 

Title card Format(20A4) 

********** 



col 1-80 20A4 Title tabular and graphics title; 

last two characters must be S' and 
immediately next to last character of 
actual title, recommend only twenty 
to thirty characters for current 
character height setting 

$$$$$$$$$$$$$$$$$$$$$$$$$$$$$S$S$$$$$$$$$$$S$$$$SSS$$$S$SSSS$$S$$S$SSSS$ 

General card 1 Format(3l5) 

********** 
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col 1-5 15 



NXT 



total x/c or y/c coordinates; 

number of x/c coordinates must equal 
number of y/c coordinates (i.e. 
they are paired) 
col 6-10 15 ITFS transition flag setter; 

set equal to 1 transition calculated 
by Michel's method which is given 
by Eq.(6.1.1) in Ref. 2 set equal 
to 0, transition is input XCTR 
(General Card 2 col 41-50) 

col 11-15 15 BC boundary conditions format; 

Set equal to 1 when inputting 
velocity ratios, UE, on General card 
2, col 31-40 , set to any integer 
value other than 1 when inputting 
the Coefficient of Pressure ,Cp, 
instead of UE 

^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ 
General card 2 Format(6F10.0) 



col 1-10 F10.0 



col 11-20 F10.0 



col 21-30 F10.0 



Col 31-40 F10.0 
Col 41-50 F10.0 



P2(1 ) 



DETA(1) 



VGP 



RC 

XCTR 



Pressure gradient parameter; 

value of m at stagnation point, set 
equal to 1.0 

initial normal step size,h; 

(i.e. is in the neta direction), 0.01 
is a good value 
Variable grid parameter, k; 

defines the ratio of successive step 
sizes in the neta direction , 1.14 is 
a good value, the higher the value 
the fewer calculations performed in 
the neta direction 
Reynolds number of chord, 
transition input; 
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x/c coordinate of transition from 
laminar to turbulent, required if 
ITFS is 0 (General card 1,Col 6-10) 
$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$s$$$s$$$$$$$$$$$$$$$$$s$$$$s$$$$$$$$ 
Data card 2 Format(3F10. 0) 

********** 



col 1-10 F10.0 XC x/c coordinate; 

starts at stagnation point, or leading 
edge, must be in increasing size 
col 11-20 F10.0 YC y/c coordinate; 

must correspond to the x/c coordi- 
nate 

col 21-30 F10.0 UE Velocity ratio; 

that is the velocity at edge of 
boundary layer divided by velocity 
at infinity, must correspond with 
location (x/c, y/c) and with BC 
assigned a value 1 (General card 1, 
col 11-15), Coefficient of Pressure 
may be inputted instead if BC is a 
value other than 1 



2. jn^uts [or CE_BCAS 



Inputs to CEBCAS are essentially the same as CEBPLT except for the 
run card. The run card is modified as follows to input other parame- 
ters needed to calculate losses and drag coefficients for cascades. 

$$$$s$$$$$$$$$$$$$$$$$$$$$$$s$$$s$$$$ss$s$$s$ssssssssss$$sss$$s$$s$$$$$$ 

Run card Format( 1 10, 3( F10. 0) ) 



col 1-10 15 N 

col 11-20 F10.0 BTA1 

col 21-30 F10.0 BTA2 

col 31-40 F10.0 SIGMA 



the number of runs (data files) 
the inlet & (degrees) 
the outlet & (degrees) 
o 
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B. OPERATION OF CEBPLT AND CEBCAS 



First you must obtain the appropriate fortran program. Link to the 
AERO disk and send appropriate program to your disk. After you got 
a copy of CEBPLT or CEBCAS on your disk you must compile the 
fortran program(s) on an appropriate fortran compiler. FILEDEF's (file 
definitions) must be defined for data input and for tabular output. An 
input data file must be created as defined. Execute DISSPLA at a 
TEK618 graphics terminal for graphics output. For specific details see 
appropriate computer manuals at the Winston Churchill Computer 
Center. 
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